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SECTION  1  TRANSFER  FUNCTION  APPROXIMATION  PROGRAM 


SCOPE 

The  Transfer  Function  Approximation  Program  (TFAP)  is  an  IBM  7090 
FORTRAN  program  which  will  compute  the  approximate  transfer  function  for  any 
linear  pressure  transducer  system  from  the  system's  time  response  to  a  step 
input.  Input  data  to  the  program  contains  information  which  specifies  test  condi¬ 
tions,  information  which  defines  the  data  range  and  density  for  the  approximation, 
and  information  which  specifies  the  amplitude  of  the  time  response  at  each  of  a 
number  of  equally  spaced  sample  points  in  the  time  domain.  Primary  output  data 
from  the  program  contains  the  amplitude  and  phase  of  the  approximate  transfer 
function  for  the  transducer  system  at  equally  spaced  points  in  the  frequency 
domain. 


With  minor  modifications,  the  program  is  adaptable  to  any  problem  which 
requires  the  computation  of  a  linear  system's  transfer  function  from  its  time 
response  to  a  known  input.  The  program  has  been  tested  with  time  response  data 
for  an  analytically  known  two  degree-of-freedom  system  excited  by  a  step  input. 
The  results  of  these  tests  indicate  that  the  accuracy  of  the  approximate  transfer 
function  is  determined  by  the  data  sampling  frequency  and  the  data  range.  For 
sampling  frequencies  greater  than  twenty  times  the  highest  frequency  of  interest 
and  a  data  range  which  includes  the  time  response  out  to  a  point  at  which  the  time 
response  has  reached  99%  of  its  final  value,  the  approximate  transfer  function 
was  found  to  be  within  .5%  of  its  known  value. 

Experimental  results  indicate  that  the  time  response  data  can  be  trun¬ 
cated  at  a  time-point  at  which  the  time  response  has  not  yet  reached  99%  of  its 
final  value  and  reasonably  accurate  results  can  still  be  expected.  If  data  is 
truncated,  the  experimental  results  indicate  that  the  maximum  error  can  be 
expected  to  occur  at  frequencies  in  the  neighborhood  of  the  system's  natural 
frequencies.  In  general,  the  maximum  error  will  occur  in  the  amplitude  charac¬ 
teristic  of  the  transfer  function  if  the  input  data  is  truncated. 

The  program  is  written  to  anticipate  data  truncation  and,  (if  the  user  so 
elects)  will  seek  a  point  for  truncation  which  satisfies  the  following  criteria: 

(1)  The  point  is  on  the  positive  going  edge  of'an  oscillation, 

(2)  The  point  is  that  point  for  which  the  amplitude  of  the  time 
response  is  near  the  average  (d-c)  value  of  the  time  response 
data, 

(3)  The  total  number  of  data  points  is  odd. 

The  first  criterion  assures  that  an  integral  number  of  cycles  of  oscillation  exist 
in  the  data  so  that  the  computed  average  value  of  the  time  response  data  is 
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based  on  an  equal  number  of  positive  and  negative  oscillations  about  the  average 
value.  The  second  criterion  Is  required  by  the  assumptions  made  in  the  mathe¬ 
matical  relations  which  define  the  approximation  procedure,  l.e. ,  the  last  data 
point  corresponds  to  the  final  value  of  the  time  response  with  all  oscillations 
dampened  out.  The  third  criterion  is  imposed  by  Simpson's  Rule  which  requires 
an  even  number  of  intervals,  i.e. ,  an  odd  number  of  data  points,  in  the  sampled 
time  response.  If  the  user  elects  to  specify  the  final  value  of  the  transient 
response,  it  is  assumed  that  the  first  two  criteria  have  been  satisfied.  The  pro¬ 
gram  will  always  satisfy  the  third  criterion  by  discarding  the  last  data  point  if  the 
given  number  of  data  points  is  an  even  number. 

The  program  removes  zero-offset  in  the  time  response  data  prior  to  the 
computation  of  the  approximate  transfer  function.  All  input  data  is  echo-checked 
by  the  program  both  before  and  after  the  final  value  has  been  determined  and  the 
zero-offset  has  been  removed. 

^Program  running  time  can  be  estimated  from  the  relation 

Running  Time  =  (NBRPTI)  (NBRPTO)  (1.7)  (lO"^)  minutes 

where  NBRPTI  is  the  number  of  time  response  data  points  used  and  NBRPTO  is  the  * 
number  of  frequency  response  data  points  to  be  computed.  For  1501  time  response 
data  points  and  201  frequency  response  points,  running  time  is  about  5  minutes. 

The  program  is  dimensioned  for  2000  input  data  points  and  an  unlimited 
number  of  output  points  and  utilizes  approximately  6000  words  of  storage  in¬ 
cluding  data  storage . 

PREPARATION  OF  INPUT  DATA 

Input  data  for  the  program  contains  three  output  identification  cards,  one 
test  data  card,  one  input/output  parameter  card,  and  the  time  response  data  cards. 
The  detailed  formats  for  these  cards  are  described  in  the  Input  Card  Format 
section  of  this  manual.  The  procedure  for  data  preparation  starts  with  the  FSADC 
transient  waveform  digitizing  operation.  Preliminary  editing  of  the  time  response 
data  is  necessary.  The  objectives  of  this  editing  are  (1)  to  remove  FSADC  output 
cards  which  correspond  to  the  oscilloscope  trace  prior  to  the  arrival  of  the  step 
pressure  function  at  the  transducer  and  (2)  to  correct  erroneous  data  points  which 
correspond  to  transparency  imperfections  or  random  errors. 

If  the  transparency  is  treated  by  applying  Kodak  Opaque  paint  to  trans¬ 
parency  imperfections  and  to  that  portion  of  the  transient  which  corresponds  to 
the  base  line,  then  the  digitizing  process  will  result  in  few  if  any  erroneous  data 
points.  It  will  still  be  necessary  to  remove  that  portion  of  the  data  which  corres¬ 
ponds  to  the  transparency  base  line,  however  all  such  data  points  will  digitize 
as  0000  for  the  amplitude.  Since  the  FSADC  prints  a  paper  tape  as  It  punches 
digitized  data  point  cards,  the  operator  can  readily  determine  the  x-coordinate  of 
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the  point  at  which  the  step  function  is  sensed.  All  data  point  cards  prior  to  that 
point  must  be  removed  from  the  card  deck  produced  by  the  FSADC.  The  total 
number  of  cards  removed  should  be  carefully  counted  so  that  the  total  number  of 
input  data  cards  to  the  program  can  be  specified.  The  FSADC  system  always 
digitizes  a  point  at  x  =  000  and  x  =  1000  ,  however  the  latter  will  appear  as 
000  on  the  corresponding  data  card  due  to  FSADC  limitations.  The  total  number 
of  data  points  will  be  101 ,  201 ,  501 ,  or  1001 ,  depending  upon  FSADC  switch 
settings  during  the  digitizing  phase. 

Assuming  that  the  transparency  was  treated  as  previously  specified,  the 
data  editing  procedure  is: 

(1)  Remove  and  count  those  data  point  cards  which  correspond  to  the 
base  line  of  the  waveform. 

(2)  Determine  the  number  of  cards  which  remain  by  subtracting  the 
number  of  cards  removed  from  the  number  of  points  digitized.  For 
example,  if  24  cards  were  removed  and  the  number  of  points 
digitized  was  501 ,  the  number  of  cards  which  remain  would  be 
477. 

(3)  Scan  the  printed  tape  for  erroneous  data  points.  The  tape  will 
contain  an  asterisk  (*)  beside  the  values  corresponding  to  such 
points  and  the  corresponding  datacard  will  be  blank  in  column  11. 
(Valid  data  cards  will  contain  a  9-punch  in  column  11 ,)  Since 
such  erroneous  points  arise  because  either  the  FSADC  sensed  two 
or  more  traces  when  it  scanned  the  x-position  or  the  FSADC  failed 
to  sense  any  trace,  the  user  must  use  judgment  in  correcting  such 
points.  Frequently,  no  correction  will  be  necessary  if  the  FSADC 
sensed  the  second  trace  at  a  point  which  was  above  the  actual 
waveform  trace.  By  examining  the  values  digitized  adjacent  to 
the  erroneous  point,  the  user  can  decide  whether  or  not  the  point 
does  require  correction.  If  correction  is  necessary,  linear  inter- 
polafton  (using  the  amplitude  values  at  adjacent  points)  should  be 
applied  in  estimating  the  correct  value.  Cards  which  require 
correction  should  be  removed  from  the  data  deck  and  marked  with 
the  corrected  value.  Preferably,  the  corrected  data  card  should 
be  keypunched  immediately  and  the  card  reinserted  in  its  proper 
position  in  the  data  deck. 

(4)  The  edited  data  deck  should  be  listed  prior  to  its  use  as  input 
data  to  the  program  and  this  listing  checked  for  missing  data 
points,  card  handling  errors,  and  keypunch  errors.  Upon  verifi¬ 
cation  of  this  deck,  the  user  should  mark  the  deck  with  such  in¬ 
formation  as  the  number  of  cards  in  the  deck  and  any  identifi¬ 
cation  which  is  considered  to  be  adequate  for  uniquely  identifying 


the  deck.  The  first  and  last  cards  of  the  deck  should  be  marked 
"1st  Card"  and  "Last  Card".  This  deck  will  be  referred  to  as  the 
"FSADC  Deck". 

Upon  completing  the  preceding  steps,  the  user  should  prepare  the  follow¬ 
ing  header  cards:  • 


(1) 

OID-1 

(2) 

OID-2 

(3) 

OID-3 

(4) 

STTD-1 

(5) 

DD-1 

The  format  and  content  of  these  header  cards  are  defined  in  the  Input  Data  Card 

Format  section.  It  is  important  that  the  user  observe  the  units  of  each  quantity 

specified  on  the  header  cards,  position  each  quantity  within  the  card  columns 

(field)  assigned  to  that  quantity,  and  keypunch  the  decimal  point  with  each 

quantity.  The  location  of  a  quantity  within  its  field  is  arbitrary  if  the  decimal 

point  is  punched.  The  only  quantity  which  must  not  have  a  decimal  point  punched 

is  NBRPTI,  i.e. ,  the  number  of  points  supplied  in  the  FSADC  Deck.  The  value 

specified  for  NBRPTI  must  be  positioned  in  columns  7-10  of  the  DD-1  card  so  that 

the  assumed  decimal  point  position  is  to  the  right  of  column  10.  ^ 

1 

INPUT  DATA  CARD  FORMAT 

Input  data  cards  to  the  Transfer  Function  Approximation  Program  are:  » 

• 

(1)  Output  Identification  Cards  —  (OID-1 ,  OID-2,  ^ID-S) 

(2)  Shock  Tube  Test  Data  Card  —  (S,TTD-1) 

(3)  Digitized  Data  Card  —  (DD-1) 

(4)  Flying  Spot  Analog-to-Digital  Converter  Cards  —  (FSADC) 

Detailed  formats  for  these  cards  are  specified  in  the  Input  Data  Card  Format  Lay¬ 
out  tables.  A  brief  description  of  the  purpose  of  each  of  the  above  card  types 
follows: 

(1)  OID-1 ,  OID-2.  OID-3  These  cards  provide  the  means  for  identi¬ 
fying  each  output  page  with  any  arbitrary  identification.  Such 
information  as  transducer  manufacturer,  transducer  serial  number, 
date,  remarks,  etc. ,  may  be  entered  on  these  cards. 

(2)  STTD-1  This  card  specifies  information  relative  to  the  shock  * 

tube  test  conditions  and  transient  waveform  recording  parameters. 

The  program  utilizes  this  information  in  computing  shock  wave 

velocity,  Mach  number,  pressure  step  size,  and  data  sampling  • 

rate. 
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(3)  DD-1  This  card  specifies  information  relative  to  the  digitizing 
phase  parameters  and  range  of  transfer  function  approximation. 

(4)  FSADC  These  cards  are  the  output  cards  from  the  FSADC 
system  that  remain  after  preliminary  editing  by  the  user.  Each 
card  represents  a  digitized  data  point  on  the  transient  waveform 
and  contains  the  transducer  system  identification  in  addition  to 
the  X  and  y  coordinates  of  the  point. 

OUTPftT  DATA  CARD  FORMAT 

Output  cards  from  the  Transfer  Function  Approximation  Program  are: 

(1)  Approximate  Transfer  Function  Header  Card  —  (ATFH) 

♦ 

(2)  Approximate  Transfer  Function  Cards  —  (ATF) 

Detailed  formats  for  these  cards  are  specified  in  the  Output  Data  Card  Format 
Layout  tables.  A  brief  description  of  the  purpose  of  each  of  the  above  card  types 
follows: 

(1)  ATFH  This  card  contains  transducer,^identification,  the  number 
of  ATF  cards,  the  bandwidth  covered  by  the  ATF  cards,  the  fre¬ 
quency  increment  between  sampled  points,  and  the  allowable  time 
interval  for  approximating  input  time  functions. 

(2)  ATF  Each  card  defines  the  approximate  transfer  function  for  the 

transducer  at  a  discrete  frequency  in  the  band  width  specified  on 
the  ATFH  card.  # 

The  above  cards  are  utilized  as  input  data  to  the  Input  Time  Function  Approxi¬ 
mation  Program.  Their  contents  also  appear  in  the  output  listings  generated  by 
the  Transfer  Function  Approximation  Program. 

PROGRAM  OPERATING  PROCEDURES 

The  order  of  input  data  cards  for  the  Transfer  Function  Approximation 
Program  is  (for  each  transducer): 

(1)  OID-1  * 

(2)  OID-2 

(3)  OID-3 

(4)  STTD-1 

(5)  DD-1 

(6)  FSADC  Deck 
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OUTPUT  DATA  CARD  FORMAT  LAYOUT 


1 


1 


**  Unlisted  columns  are  blank 


The  above  card  deck  is  referred  to  as  a  Transducer  Data  Deck.  An  unlimited 
number  of  these  decks  can  be  processed  with  a  single  run  on  the  IBM  7090  by 
stacking  the  decks  in  the  desired  order  of  processing.  The  fiist  deck  must 
always  be  preceded  by  a  card  with  an  asterisk  (*)  punched  in  column  1  and  the 
word,  DATA,  punched  in  columns  7-10.  The  last  deck  must  be  followed  by  an  end- 
of-file  card,  i.e. ,  a  card  with  a  7  and  8  punch  in  column  1  and  the  end-of-file 
card  must  be  followed  by  a  card  with  an  asterisk  in  column  1  and  the  words,  END 
TAPE,  in  columns  7-14. 

The  user  should  consult  with  the  computing  installation  staff  for  the  card 
formats  they  have  established  for  the  FMS  (Fortran  Monitor  System)  date  and 
identification  cards.  These  two  cards  must  precede  the  TFAP  Binary  Deck  (i.e. , 
the  TFAP  program  deck  in  column-binary  form).  The  user  must  assemble  the  pre¬ 
ceding  cards  and  decks  into  the  following  order  for  a  computer  run: 

(1)  FMS  Date  Card 

(2)  FMS  Identification  Card 

(3)  TFAP  Binary  Deck 

(4)  *  Data  Card 

(5)  Transducer  No  1  Data  Deck 

(6)  Transducer  No  2  Data  Deck 

(7)  Transducer  No  N  Data  Deck  « 

(8)  End-of-File  Card 

* 

(9)  *  End  Tape  Card 

The  TFAP  program  reads  all  input  from  logical  tape  unit  No  5,  writes  all 
output  for  listing  purposes  on  logical  tape  unit  No  6,  and  writes  all  output  for 
card  punching  on  logical  tape  unit  No  11 .  All  input  and  output  is  in  BCD  form. 

The  program  will  write  an  End  File  mark  behind  each  output  on  logical  tape  unit 
No  11.  The  total  number  of  such  End  File  marks  is  equal  to  the  total  number  of 
Transducer  Data  Decks  submitted  by  the  user  for  the  computer  run.  The  user  must 
specify  the  preceding  information  to  the  computing  installation  staff. 

The  TFAP  program  makes  use  of  the  on-line  printer  during  processing. 
When  called  upon&to  compute  the  estimated  final  value  of  the  transducer  time 
response  data,  the  program  will  print  the  results  of  each  value  estimated  during 
the  iteration  process  until  the  estimation  criteria  has  been  satisfied.  Generally, 
no  more  than  15  to  20  lines  of  output  should  occur  during  the  estimating  phase, 
however  the  user's  experience  with  the  program  will  determine  what  is  reasonable. 
The  Program  Organization  section  should  be  consulted  for  further  details  on  the 
estimating  process.  Following  a  successful  estimation,  the  program  will  print 
on-line  an  estimated  running  time  figure  for  the  Transducer  Data  Deck  currently 
being  processed.  The  u^er  should  advise  the  computing  installation  staff  of  the 
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approximate  value  to  expect  for  each  transducer  since  an  error  in  card-handling 
or  data  preparation  will  usually  result  in  either  a  ridiculous  time  estimate  or 
excessive  processing  time  by  the  program.  This  information  will  tend  to  reduce 
the  possibility  of  unproductive  computer  runs  due  to  card-handling  and/or  key¬ 
punching  errors. 

Upon  completion  of  a  transfer  function  approximation  for  each  Trans¬ 
ducer  Data  Deck,  the  program  will  print  on-line  that  it  has  completed  the 
approximation.  This^rint  should  occur  within  the  time  period  which  began  frith 
the  printing  of  the  estimated  running  time.  ^  ^ 

PROGRAM  ORGANIZATION 

It  is  assumed  that  the  user  is  familiar  with  the  Mathematical  Analysis 
section  of  this  report.  The  basis  for  the  computational  procedures  utilized  in  the 
Transfer  Function  Approximation  Program  are  developed  and  discussed  in  that 
section  in  Volume  I. 

The  block  diagram  for  the  program  is  presented  in  this  section.  The 
reader  should  refer  to  the  Program  Listing  section  for  the  program  details  which 
represent  each  block.  The  following  conventions  are  established  as  an  aid  to 
associating  the  block  diagram  with  the  corresponding  FORTRAN  statements: 

(1)  Where  a  block  is  numbered,  that  number  corresponds  to  the 
FORTRAN  statement  number  which  begins  the  sequence  of  oper¬ 
ations  defined  in  the  block.  The  block  number  will  Be  oiftside 
the  block  at  the  upper  left  hand  corner. 

(2)  Input  and  output  blocks  will  contain  the  corresponding  FORTRAN 
format  statement  number  in  the  lower  left  hand  corner  of  the  block. 
The  number  of  the  logical  tape  unit  bei|g  used  as  the  input  or  out¬ 
put  unit  will  appear  in  the  lower  right  Ipnd  corner  of  the  block. 

(3)  •  Each  FORTRAN  statement  card  is  numbered  in  columns  75  through 

80.  The  beginning  and  ending  card  number  for  %  sequence  of 
operations  represented  in  the  block  diagram  will  appear  outside 
the  corresponding  block  at  the  upper  and  lower  right  hand  corners 
of  the  block,  respectively.  * 

The  major  computational  processes  in  the  program  are: 

(1)  Final  Value  Estimation 

(2)  Transfer  Function  Approximation 

The  program  always  checks  to  see  if  the  number  of  raw  time  response  data  points 
is  an  odd  number.  If  not,  the  last  data  point  is  discarded  and  the  total  input  data 
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point  count  is  reduced  by  one.  The  program  then  checks  the  value  given  by  the 
user  on  the  DD-1  card  for  the  final  value.  If  the  given  value  is  zero,  the  program 
assumes  that  the  final  value  must  be  estimated.  Prior  to  beginning  the  estimation, 
the  program  checks  for  a  positive  slope  on  the  portion  of  the  time  response  defined 
by  the  last  two  time  response  po^ts.  If  the  slope  is  not  positive,  the  program 
discards  the  last  two  data  points  and  reduces  the  total  input  data  point  count  by 
two.  The  slope-check  is  repeated  and  points  are  discarded  by  two's  until  the  pro¬ 
gram  determines  that  the  last  two  points  do  define  a  positive  slope.  The  program 
then  determines  the  average  value  of  the  remaining  raw  time  response  data  by 
numerically  integrating  the  data  and  dividing  by  the  time  range  defined  by  the 
data.  The  average  value  is  compared  to  the  last  time  response  data  point.  If  the 
average  value  differs  from  the  iast  time  response  data  point  by  an  amount  which 
is  less  than  the  difference  between  the  last  data  point  and  the  second  from  the 
last  data  point,  the  program  assumes  that  the  final  value  of  the  time  response 
data  is  the  computed  average  value.  If  this  difference  criteria  is  not  satisfied, 
the  program  discards  the  last  two  data  points  and  returns  to  the  section  which 
checks  for  positive  slope.  The  program  will  print  on-line  (1)  the  value  of  the 
second  from  last  data  point,  (2)  the  computed  average  value,  and  (3)  the  value  of 
the  last  data  point,  prior  to  returning  to  the  slope-check  section.  Since  all 
numerical  Integration  performed  in  the  program  is  based  on  Simpson's  Rule,  data 
points  must  be  discarded  by  two's  to  preserve  the  "even  number  of  intervals" 
requirement  (i.e. ,  an  odd  number  of  points  are  required  when  Simpson's  Rule  is 
used) .  This  method  of  estimating  the  final  value  resulted  in  an  estimated  final 
value  of  .9991  versus  an  actual  final  value  of  .9998  or  about  .0  7%  error  for  the 
known  two-degree-of-freedom  test  case. 

The  section  of  the  program  which  computes  the  approximate  tranter 
function  is  based  on  the  relations  defined  in  Equations  (44)  through  (50)  of  the 
Mathematical  Analysis  section  of  the  report.  The  integrals  are  evaluated  numeri¬ 
cally  by  applying  Simpson's  Rule.  It  should  be  noted  that  the  phase  angle  for  the 
approximate  transfer  function  is  determined  in  the  program  from  the  relation 

-1  * 

<t>  =  tan  Y/X  -  270° 

* 

and  not  the  relation  defined  in  Equation  (^).  Experimental  results  with  the  two 
degree-of-freedom  system  indicate  that  a  phase  error  will  occur  if  the  0/2  angle 
is  included,  and  that  this  phase  error  is  approximately  equal  to  the  angle,  G/2. 

It  appears  that  the  source  of  this  error  is  in  the  approximation  made  for  the 
Fourier  Transform  of  a  step  function.  The  program  also  includes  ^provision  for 
completing  the  angle,  tan"^Y/X  ,  when  it  exceeds  270°.  This  provision  forces 
the  phase  angle  to  always  be  developed  as  a  negative  angle. 

A  major  part  of  the  program  is  concerned  with  the  details  of  the  house¬ 
keeping  that  are  necessary  for  obtaining  the  output  formats  and  data.  It  should 
be  noted  that  the  program  echo-checks  the  raw  time  response  data  and  follows 
this  with  the  processed  time  response  data,  i.e. ,  the  time  response  data  with 
zero  offset  removed. 


TRANSFER  FUNCTION  APPROXIMATION  PROGRAM  BLOCK  DIAGRAM 
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The  preparation  and  editing  of  time  response  data  Is  a  function  of  the 
system  being  used  to  generate  and  digitize  this  data.  If  the  FSADC  is  being 
used,  the  preparation  and  editing  procedures  are  the  same  as  those  described  in 
the  TFAP  write-up.  Regardless  of  the  system  being  used  to  obtain  this  data,  the 
card  format  for  time  response  data  points  must  be  identical  to  the  FSADC  card 
format  which  is  also  described  in  the  TFAP  write-up. 

The  user  must  prepare  the  three  output  identification  cards  and  the 
input/output  parameter  card  in  accordance  with  the  formats  defined  in  the  Input 
Data  Card  Format  section. 


INPUT  DATA  CARD  FORMAT 


Input  data  cards  to  ITFAP  are: 

(1)  Output  Identification  Cards  —  (OID-1,  OID-2,  OID-3) 

(2)  Approximate  Transfer  Function  Header  Card  —  (ATFH) 

(3)  Approximate  Transfer  Function  Cards  —  (ATF) 

(4)  Input/Output  Parameter  Card  —  (lOP) 

(5)  Time  Response  Data  Cards  —  (TRD) 


The  detailed  card  format  for  the  OID  cards  is  described  in  the  TFAP  write-up  as 
are  the  ATFH  and  ATF  cards.  The  TRB  cards  must  have  the  same  format  as  the 
FSADC  cards  described  in  the  TFAP  write-up,  however  the  contents  of  columns 
1-10  of  these  cards  should  contain  time  response  data  run  identification.  The 
lOP  card's  detailed  format  is  specified  in  the  Input  Data  Card  Layout  table.  A 
bri^  description  of  the  information  contained  on  this  card  follows.  It  specifies 
the  time  response  data  run  identification,  the  number  df  time  response  data 
points,  the  highest  frequency  component  assumed  to  exist  in  the  input  time 
function,  the  time  response  data  sampling  interval,  the  time  increment  to  be 
used  between  successive  input  time  function  points,  the  time  interval  for 
approximation  of  the  input  time  function,  and  the  conversion  constant  for  con¬ 
verting  time  response  amplitude  values  to  pressure. 


PROGRAM  OPERATING  PROCEDURES 

The  order  of  input  data  for  ITFAP  is: 


(1) 

OID-1 

(2) 

OID-2 

(3) 

OID-3 

(4) 

ATFH 
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INPUT  DATA  CARD  FORMAT  LAYOUT 
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(5) 

ATF  Deck 

(6) 

lOP 

(7) 

TRD  Deck 

Any  number  of  input  data  decks  can  be  processed  with  a  single  run  on  the  IBM 
7090  by  assembling  each  deck  in  the  above  order  and  stacking  the  assembled 
decks  in  the  desired  order  of  processing.  The  first  such^eck  must  be  preceded 
by  an  *Data  card  and  the  last  such  deck  must  be  followed  by  an  End-of-Flle 
and  End  Tape  card.  These  cards  are  described  in  the  TFAP  write-up. 

The  user  should  consult  with  the  computing  installation  staff  for  the 
card  formats  they  have  established  for  the  FMS  (FORTRAN  Monitor  System)  date 
and  identification  cards.  These  two  card^must  precede  the  ITFAP  Binary  Deck, 
i.e. ,  the  ITFAP  program  in  column-binary  form .  The  user  must  assemble  the 
preceding  cards  and  decks  into  the  following  order  for  a  computer  run; 


(1) 

FMS  Date  Card 

(2) 

FMS  Identification  Card 

(3) 

ITFAP  Binary  Deck 

(4) 

*  Data  Card 

•  •  • 

(5) 

Data  Decks 

(6) 

End-of-File  Card 

1 

(7) 

*  End  Tape  Card 

•  The  ITFAP  program  reads  all  input  data  from  ^gical  tape  unit  No  5  and  4 

writes  all  output  for  listing  purposes  on  logical  tape  unit  No  6.  All  input  and 
output  is  in  BCD  form. 

•  • 

ITFAP  makes  use  of  the  on-line  printer  during  processing.  The  pro¬ 
gram  prints  an  estimate  of  the  running  time  for  each  data  deck  prior  to  process¬ 
ing  the  data  deck.  The  user  should  advise  the  computing  installation  staff  of 
the  approximate  value  to  expect  for  each  data  deck  since  an  error  in  card  hand¬ 
ling  or  data  preparation  will  usually  result  in  either  a  ridiculous  time  estimate 
or  excessive  processing  time.  The  program  prints  an  on-line  estimate  as  pro¬ 
cessing  of  each  data  deck  is  completed  and  this  comment  should  appear  within 
the  time  which  has  lapsed  since  the  running  time  estimate  was  printed. 

PROGRAM  ORGANIZATION 

It  is  assumed  that  the  user  is  familiar  with  the  Mathematical  Analysis 
section  of  this  report  and  the  TFAP  write-up.  The  basis  for  the  computational 
procedures  utilized  in  the  approximation  of  input  time  functions  are  developed 
and  discussed  in  the  former  and  the  block  diagram  conventions  are  established 
in  the  latter. 
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The  major  computational  processes  in  the  program  are: 


(1)  Input  Function  Spectrum  Computation 

(2)  Input  Time  Function  Approximation 

The  program  always  checks  to  see  if  the  number  of  time  response  data  points  is 
an  odd  number  and  will  discard  the  last  data  point  if  not.  This  requirement  is 
due  to  the  method  of  numerical  integration  utilized  in  the  program  (Simpson's 
Rule) .  After  echo-checking  the  input  data,  the  program  removes  the  x  and  y 
zero  offset  in  the  time  response  data  and  estimates  the  running  time.  The  latter 
will  appear  in  an  on-line  comment.  This  estimate  is  conservative  and  if  the 
running  time  exceeds  this  figure,  it  is  reasonable  to  assume  that  a  malfunction 
has  occurred  due  to  incorrect  data  preparation. 

The  program  computes  the  frequency  spectrum  for  the  input  time  function 
by  computing  the  frequency  spectrum  for  the  time  response  data  over  the  range 
specified  by  the  user,  and  multiplying  this  spectrum  by  the  inverse  transfer 
function.  Since  th^,  transfer  function  is  known  at  discrete  points  only,  the  time 
response  spectrum  is  computed  at  the  same  points  only.  Usually,  the  magnitude 
of  the  high  frequency  components  of  the  ir>put  time  function  will  be  small  relative 
to  the  lower  frequency  components.  If  so,  computation  time  can  be  reduced  By 
the  user  through  specifying  the  limiting  frequency  at  the  lowest  feasible  value. 
The  user  may  view  this  part  of  the  program  as  a  low-pass  filter  whose  band¬ 
width  is  equal  to  the  difference  between  FREQLM  (Specified  on  lOP  header  card 
in  ITFAP  data)  and  BEFREQ  (originally  specified  on  the  DD-1  header  card  in 
TFAP  data) .  It  is  interesting  to  note  that  the  program  can  be  readily  modified  to 
act  as  a  high-pass,  band-pass,  or  notched  filter  with  sharp  roll-off  character¬ 
istics.  The  program  evaluates  the  frequency  spectrum  for  the  time  response 
data  by  numerical  computation  of  the  values  defined  in  Equations  (44)  and  (45) , 
over  the  bandwidth  specified  in  the  input  data  and  in  increments  defined  by  the 
approximate  transfer  function  data.  As  each  value  of  X(juj)  is  computed,  the 
program  divides  its  amplitude  by  the  amplitude  of  the  approximate  transfer 
function  at  that  frequency  to  obtain  the  amplitude  of  F(ja)) ,  The  phase  of  F(jw) 
is  computed  by  computing  the  difference  between  the  phase  of  X(jw)  and  the 
phase  of  the  approximate  transfer  function.  These  two  results  are  used  to  re¬ 
solve  F(jw)  into  its  real  and  imaginary  parts.  The  program  stores  these  two 
parts  for  use  during  the  input  time  function  approximation  phase. 

Upon  completing  the  computation  of  the  spectrum  of  F(jcj)  ,  the  pro¬ 
gram  proceeds  to  the  section  which  computes  the  approximate  input  time  begin¬ 
ning  at  t  =  0  ,  The  computation  is  based  on  the  relations  given  in  Equations 
(61)  and  (62).  Both  equations  are  evaluated  by  numerical  integration  and  the 
results  from  the  two  equations  will  be  approximately  equivalent  if  F(jw)  is 
Fourier  transformable.  In  order  to  aid  in  the  interpretation  of  results,  each 
equation's  result  will  appear  in  the  output  listing  as  will  the  nearest  value  of 
time  response.  Equation  (61)  is  referred  to  as  the  "Cosine  Integral"  and 
Equation  (62)  is  referred  to  as  the  "Sine  Integral". 


The  user  should  note  that  the  zero  offset  is  removed  from  the  time 
response  data  prior  to  its  use  In  the  computation,  hence  the  conversion  con- 
,  stant  which  Is  specified  by  the  user  should  be  given  on  a  pressure  per  unit 

amplitude  basis.  • 

t  « 
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INPUT  TIME  FUNCTION  APPROXIMATION  PROGRAM 
BLOCK  DIAGRAM 
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HERCO  PROJECT  P-106-9-  2  AHM  -  2-21-63  TFAP-02  0001 

XEQ  0002 

FORTRAN  0003 

TRANSFER  FUNCTION  APPROXIMATION  PROGRAM  0004 

PROGRAMMED  FOR  0005 

6593D  TEST  GROUP ( DEVELOPMENT )  -  AIR  FORCE  FLIGHT  TEST  CENTER  -  0006 

EDWARDS  AIR  FORCE  BASE*  CALIFORNIA  0007 

PROGRAMMED  BY  0008 

HOUSTON  ENGINEERING  RESEARCH  CORPORATION  0009 

UNDER  0010 

CONTRACT  NUMBER  AF  04(611)  8199  0011 

DIMENSION  XT(2000)«  HEDI NO ( 36 ) »XCORD < 2000 ) • I DOUT ( 2 ) (XARRAY ( 300 ) »  0012 

lYARRAYOOO)  0013 

1  READ  INPUT  TAPE  5 . 10 1 . ( HEDI NG ( I ) » I »1 ♦ 36 )  0014 

101  FORMAT! 12A6/12A6/12A6)  0015 

20  READ  INPUT  TAPE  5 . 102 .P 1 f PA . TEMP .OIST ♦ TI M I  NT .GAMMA »R » VERSEN .HORSWP  0016 

102  FORMAT! 9F8.0)  0017 

30  READ  INPUT  TAPE  5 .103.NBRPTI tXTFIN.  BEFREQ.FIFREQ.FREQIN  0018 

103  FORMAT!  no, 4F10.0)  *  0019 

READ  INPUT  TAPE  5 , 1 1 3 , I DOUT ! 1 ) . I DOUT ! 2 ) .XCORD 1  A .X T!  1 )  0020 

113  FORMAT!2I4,9X,F3.0,F10.0)  0021 

•  40  READ  INPUT  TAPE  5 » 104. ( XCORD ( I ) ,XT ! I > , I *2 .NBRPT I )  •  0022 

104  FORMAT! 15X,F5.0.F10.0)  0023 

DELX  =  XCORD(2)  -  XCORD! D  0024 

C  VELOCITY  -  MACH  NBR  -  PRESSURE  STEP  COMPUTATION  0025 

r  VEL  =  DIST*1.E6/TIMINT  0026 

FMACH  =  VEL/SQRTF! 32.2»GAMMA»R*!TEMP+460. ) )  0027 

FM2  =  FMACH*FMACH  0028 

PSTEP  =  (Pl+PA)*(2.3333333»!FM2-l.)»(4**FM2+2.n/(FM2+5.)  0029 

IBASE  =  0  0030 

23  IFdBASE  -  NBRPTI)  24,28.28  *  0031 

24  NPAGE  =  NBRPTI  -  IBASE  0032 

IF(NPAGE  -  300)  25,26,26  0033 

25  LIM2  =  NPAGE  0034 

GO  TO  27  0035 

26  LIM2  =  300  0036 

27  CALL  HEDOUT!  IBASE, HEOING)  «  0037 

DO  188  I  =  1,  LIM2  0038 

1C  =  I  +  IBASE  0039 

XARRAY!!)  =  XCORD!K)  0040 

188  YARRAY(I)  »  XT!  PC)  0041 

WRITE  OUTPUT  TAPE  6 , 1 16 , I OOUT ! 1 ) , IDOUT ! 2 ) , ( XARRAY ! I ) , YARRAY ( I ) , I -1  0042 

1,LIM2)  0043 

116  FORMAT! 1H0.49HRAW  INPUT  TIME  RESPONSE  DATA  FOR  TRANSDUCER  I*0.  ,  0044 

12I4,/6!2F10,0) )  0045 

IBASE  *  IBASE  +  LIM2  0046 

GO  TO  23  0047 

28  CALL  HEDOUT! l.HED I NG)  0048' 

«  C  FINAL  VALUE  ESTIMATION  0049; 


CXTFIN«XTFIN 

NBRPT-NBRPTI 

IF(NBRPT  -  (NBRPT  /2 ) *2 » 301 t300#301 

300  NBRPT  ■  NBRPT  -1 

301  CONTINUE 
IF(XTFIN)50,200.50 

200  PRINT  208f rDOUT(l>.IDOUT(2) 

208  FORMAT (1H1,42HC0MPUTED  FINAL  VALUES  FOR  TRANSDUCER  I .D*  .2I4/1H0# 
llOX.lOHXTt LAST-2) . 7X . ISHAVERASE  VALUE » 12X.BHXT ( LAST ) ) 

202  BASE«NBRPT-1 
LIM=BASE 

DIF1=XT( NBRPT) -XT (LIM-1) 

IF(DIF1)203.201.201 

203  NBRPT=NBRPT-2  • 

GO  TO  202 

201  SUMOD-0, 

SUMEV-0, 

204  DO  205  I»2.  LIM.  2 
SUMODiSUMOD+XT( I ) 

205  SUMEV=SUMEV+XT( I+l) 

SUMEV=SUMEV-XT( NBRPT ) 

CXTFINal XT( 1 )+XT(NBRPT )+4.*SUMOI>+2«*SUMEV)/(3.*BASE) 
DIF2=ABSF(XT(NBRPT )  -  CXTFIN) 

IF(DlPl-DIF2)206»  50.  50 

206  NBRPT=NBRPT-2 

PRINT  207. XTtNBRPT) .CXTFIN»XT(NBRPT+2)  • 

207  FORMAT! IH  .  3F20.5) 

GO  TO  202 

50  WRITE  OUTPUT  TAPE  6.  444 
444  FORMAT! IHO/lHO/lHO/lHO/lHO/lHO/) 

ECHO  CHECK  AND  PRELIMINARY  OUTPUT 

WRITE  OUTPUT  TAPE  6 . 106 .P 1 .P A .TEMP .DI ST .T IMI NT .R .GAMMA . VERSEN .HORS 
IWP 

J06  FORMAT! lH0.48X.23HeCH0  CHECK  OF  TEST  DATA/ IHO . 23HTEST  SECTION  PRES 
ISURE  =.F10.3<5H  PSIG.21X.22HATMOSPHERIC  PRESSURE  =.F10.3i5H  PSIA/1 
2H0.13HTEMPERATURE  =.F10.3.6H  DEG-F .30X . 33HD I  STANCE  BETWEEN  VELOCIT 
3Y  GAGES  =.-F10.3.5H  FEET/IHO.  15HTIME  INTERVAL  =.F10.3.13H  MICROSECO 
4NDS.21X.23HSPEC1FIC  GAS  CONSTANT  =.F10.3.17H  FT-LBS/LBM-DEG-R/IHO. 
521HSPECIFIC  HEAT  RATIO  =.F10i3.27X»23H  VERTICAL  SENSITIVITY  «.F6.Q 
6.14H  MILL1VOLTS/CM/1H0.18HH0RI2ONTAL  SWEEP  «»F10.0.16H  MICROSECOND 
7S/CM) 

WRITE  OUTPUT  TAPE  6.115 
115  FORMAT! 1H0///50X.18HC0MPUTED  TEST  DATA) 

CONVER  =  5. E2/(CXTFIN-XT! 1 ) )*PSTEP /VERSEN  P 

WRITE  OUTPUT  TAPE  6 . 107 .VEL .FMACH . PSTEP. CONVER 
107  FORMAT!  IHO. 21HSHOCK  WAVE  VELOCITY  *.F10.3.7H  FT/SEC. 21X..13HMACH  NU 
IMBER  *.F10.3./IH0.25HREFLECTE0  PRESSURE  STEP  -.FIQ.S.AH  PSI.20X. 
224HTRANSDUCER  SENSITIVITY  ».F10.3.14H  PSI/MILLI VOLT/1H0///5OX. 
323HINPUT/OUTPUT  PARAMETERS) 

• 


0050 
0051 
0052 
0053 
0054 
0055 
0056  r'. 
0057 
0058 
00  59 
0060 
0061 
006  2 
006«  ' 
0064 
0065  ,  ■ 
0066 
0067 
0068  ,~ 
0069  " 
0070 
0071  '■ 
0072 
0073 
0074  ' 
0075  '  ' 
0076 
0077  ■' 
0077 ;i 
0077.2 
0078 
0079  ' 
008  0 
0081  ' 
0082 
0083 
0084  • 
0085  ■ 
0086 
0087  ' 
0038 
0089 
0090 
0091  ' 
0092 
0093 
0094 
0095 
0096 
0097 


DELTI  -  DELX»H0RSWP*l.E-8  0098.' 

NBRPTO-(  FlFREQ-BEFREO)/FREQm  +1*  •  0099 

IF(NBRPTO-(NBRPTO/2)*2»401»40©»401  0100 

400  NBRPTO»NBRPTO+1  .  0101." 

401  CONTINUE  0102 

FFN  ■  NBRPTO  -  1  0103 

FIFREQ  ■FFN*FREQIN  +  BEFREQ  0104. 

WRITE  OUTPUT  TAPE  6ilQ8*MBRPTI»M8RPT»DELTI  . BEFREQ .F I  FREQ. FREQ  I N  0105 

108  FORMAT (1H0,38HNUMBER  OF  INPUT  DATA  POINTS  SUPPLIED  » » 1 4 . 18X. 34HNUM  0106 

IBER  OF  INPUT  DATA  POINTS  USED  I4./1H0.24HDATA  SAMPLING  INTERVAL  0107  ' 
2-.Fll#8.13H  SECONDS  AND  .  25HTRANSFER  FUNCTION  RANG?  ■.F6*0*2H  -»F  0108  ' 

38.0.4H  CPS,12H  IN  STEPS  0F.F8.0.4H  CPS)  0109 

SAMP1=.5/DELT1  0110,' 

SAMP2».5/FREQIN  0111" 

WRITE  OUTPUT  TAPE  6 . 109 .SAMP  1 .SAMP2  0112 

109  FORMAT! 1H0,34HSAMPLING  THEOREM  FREQUENCY  LIMIT  »*F10.0.8H  CPS  AND.  0113 

138H  INPUT  TIME  FUNCTION  INTERVAL  LIMIT  =  .F10.8.  8H  SECONDS)  0114 

WRITE  OUTPUT  TAPE  6.  117 .  XT ( NBRPT  I ) »CXTF IN . XCORD ( 1 )  .XT ( 1 )  0115 

117  FORMAT! 1 HO, 22HRAW  DATA  FINAL  VALUE  »»F10.0.23H  AND  FINAL  VALUE  USE  0116, 7 
10  =,  F10,0.3X,15HZERO-OFFS£T  X  *»F10.0.1X. 15HZERO-OFFSET  Y  *  0117 

2F10.0)  0118 

NBRPTI=NBRPT  ,  0119 

«  XCBASE  «  XCORD(l)  0120 

XTBASE  -  XTU  )  0121 

CXTFIN  -  CXTFIN  -  XTBASE  0122  ' 

DO  2  I=1.NBRPTI  0123 

XCORD! I)  »  XCORD! I)  -  XCBASE  0124 

2  XT!I)  »  XT!I)  '-  XTBASE  0125  ''. 

RANGE  =  l.E8/(XCORDIN8RPTl)*HORSWP)  0125.1 

WRITE  OUTPUT  TAPE  6, 119, RANGE  0125.2 

119  FORMAT! 1H0.65HTHE  TRANSFER  FUNCTION  MAY  BE  INACCURATE  FOR  FREQUENC  0125v' 
llES  LESS  THAN,F10.0.40H  CPS  IF  TIME  RESPONSE  DATA  WAS  TRUNCATED)  0125.4 
IBASE  *0  . 

63  IF! IBASE  -  NBRPTl ) 64 ,68 ,68  0127  ' 

64  NPAGE  =  NBRPT I  “  IBASE  0128 

1F!NPAGE-300)  65,66,66  0129 

65  LIM2  *  NPAGE  0130  " 

GO  TO  67  0131 

66  L1M2  “  300  0132 

67  CALf  HEDOUT! l.HEDING)  0133, • 

DO  199  I  =  l,LrM2  0134 

K  =  I  +  IBASE  0135 

XARRAVII)  -  XCORD IK)  0136.' 

199  YARRAYIl)  »  XT(K)  0137' 

WRITE  OUTPUT  TAPE  6 , 54 , ( XARRAY ( I ) .YARRAYIl » » I ■! rUI M2 )  0138 

54  FORMAT! 1H0,78HINPUT  DATA  USED  IN  APPROXIMATION  OF  TRANSFFR  FUNCT 10  0139  ‘ 

IN  WITH  T»0  OFFSET  REMOVED  /6(2F10iG))  0140 

IBASE  *  IBASE  +  LIM2  0141  J 

60  TO  63  0142  'I 

68  ESTIME  NBRPTOiNIBRPTI  0143  | 
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o 


ESTIME  »  ESTIME  *  1.6E-5 

PRINT  151»  ESTIME»IDOUT<1)»IDOUT{2I 

151  FORMAT! 1H1»  47HTHE  NEXT  ON  LINE  P«INT-OUT  SHOULD  OCCUR  WITHIN  .F4* 

10.  24H  MINUTES  FOR  TRANSDUCER  .2I4»/1H1) 

WRITE  OUTPUT  TAPE  11.  132 . lOOUT ( 1  r*100UTT2 ) »NBRPTO»ftEFREQ»FiFREOf 
1  FREQIN.SAMP2 

152  FORMAT(2I4.I12»3F10.0»F15*8) 

TRANSFER  FUNCTION  APPROXIMATION 
INIT  =50 

FREQ  =  BEFREQ 

DO  49  J=1.NBRPT0 

ANGLE  =  6.2831S53  •FRE0*DELT1 

SUMOl  *  0, 

SUM02  =  0. 

SUMEl  »  0. 

SUME2  =0.  • 

LIH  >■  NBRPTl  -  1 
DO  3  I  -  2.LIM.2 
F'l  =  1-1 

TANGl  »  FI  *  ANGLE 

SUMOl  =  SUMOl  +  XT! I )*C0SF(TANG1 » 

SUM02  =  SUM02  +  XT ! I ) *S I NF  !T ANGl ) 

12  FI  ■=  I 

TANGl  =  FI*ANGLE 

SUMEl  =  SUMEl  +  XT! I+l )*C0SF!TANG1 » 

3  SUME2  «  SUME2  +  XT! I+l )*SINF!TANG1 ) 

SUMEl  =  SUMEl  -  XTINBRPTI )*COSFtTANGl) 

SUME2  =  SUME2  -  XT ! NBRPTl )*S INF ( TANGl ) 

10  FN  »  NBRPTl  -  1 
TANGl  *  FN»ANGLE 

YOPYNl  «  XT!1)  +  XT ! NBRPTl )*COSF( TANGl) 

Y0PYN2  *  XT! NBRPTl )*SINF! TANGl) 

SUMl  =  .33333333+!Y0PYNl  +  4.+SUM01  +  2.*SUMEl) 

SUM2  =  .33333333*! Y0PYN2  +  4.*SUM02  g  2.*SUME2) 

TANGl  =  ANGLE*! FN+.5) 

FACT=1./S1NF( .5*ANGLE) 

X  =  SUMl  -  .5*  CXTFIN  *S INF ! TANGl ) *F ACT 
Y  *  -SUM2  -  .5*  CXTFIN  *COSF TTANGl )*FACT 
AMP  =  SQRTF  !X*X+Y*Y) 

AMPLIT  =  2./  CXTFIN  *S I NF ! ANGLE*.5 )*AMP 
PHASE  =  ARCTAN!Y.X ) 

IFIPHASE  -  4.7123889)  16.16.15 

15  PHASE  =  PHASE  -  6.2831853 

16  PHASER  *  PHASE  -  4.7123889 
PHASED  =  PHASER  #  57.29578 
IF! 1NIT-50)502.501.501 

501  CALL  HED0UT!1*HED1NG)  • 

WRITE  OUTPUT  TAPE  6.110.PSTEP 

110  FORMAT! 1H0.22HUN1T  PRESSURE  STEP  IS  »F4.0»4H  PSI/IH  .ISHFREQUENCY- 
lCPSt5X«9)iANPLITU0E»5X*9HPHASE-&EGt5X«9HPHASE‘‘RA0»5X«9HTHETA-0£G»9X 


0144  0 
0145 
0146 
0147  0 
0148  " 
0149 
015  0 
0151 
0152 
0153  . 


0154 


0155 
0156  O 
0167  ' 
0158 
0159 
0160 
0161 
0162 
0163  '  ■ 
0164 
0165  0 
0166 
0167 
0168  0; 
0169 
0170 
0171  ,■> 
0172 
0173 
01740 
0175 
0176 


0177 

0178 

0179 

0160 

0161 

0182 

0163 

0184 

0185 

0186 

0187 


v-y 


0188 
0189  '  . 
0190  " 
0191 
0192  .0 
0193 


<.j 


o 
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2  »11HX-C0MP0NENT*9X.11HY-C0MP0NENT/1H  ) 

INIT*0 

502  THETA=ANGLE*57. 29578  » 

WRITE  OUTPUT  TAPE6 , 1 1 1 . FR EQ  *AMPLI  T».PHASED.  PHASER »THETA,X  tY 

111  FORMATdH  .F13.0.4F14.3.2E20*3> 

WRITE  OUTPUT  TAPEll . 1 12  *  I DOUTI 1 ) . lOOUT ( 2 ) »FREO»AMW.I T ♦PHASeO»THETA 

112  FORMAT(2r4.F12.0.3F15.5) 

I  NIT  =  INIT  +  1 

49  FREQ  =  FREQ  +  FREQIN 

PRINT  114.ID0UT(1)»100UT(2) 

114  FORMAT (34H0TRANSFER  FUNCTION  FOR  TRANSDUCER  »2I4«19H  HAS  REEN  COMP 
lUTED.) 

END  FILE  11 
GO  TO  1 
END 

FORTRAN 

SUBROUTINE  HEDOUTC  HtHEOlNG) 

DIMENSION  HEDING(36)  • 

IF(  M  )  802t801«802 

801  IPAGE  =  M 

802  IPAGE  =  IPAGE  +  1 

WRITE  OUTPUT' TAPE  6*  803 • IPAGE * ( HEOING ( 1) • I > 1 >36 ) 

803  FORMAT! IHl «57X> 5HPAGE  d2>/24X*12A6/24X>I2A6/24X>12A6) 

RETURN 

END 

FORTRAN 

FUNCTION  ARCTAN(Y*X) 

IF  (X)  1«2.3 

1  IF(Y)11«12>13 

11  R~  3.1415927 
GO  TO  5 

12  R  =  3.1415927 
GO  TO  4 

13  R  *  -3.1415927 
GO  TO  5 

2  IF(Y)  21.22.23 

21  R  »  4.7123889’ 

GO  TO  4 

22  R  »  0.  • 

GO  TO  4 

23  R  “  1.5707963 

GO  TO  4  * 

3  IF  (Y)  31.32.33 

31  R  =  -6.2831854 
GO  TO  5 

32  R  =  0.  • 

GO  TO  4 

33  R  *  0. 

GO  TO  5 

4  ARC TAN  «  R 


0194  ■■ 
0196  ' 
0196 
0197  (' 
0198  ' 
0199 
0200  (" 
0201  " 
0202 
020  3  ,  - 
0204  ' 
0205 
0206  ;  ■ 
0207  ' 
0208 
0208»:* 
0209  ' 
0210 
0211  {.■•■ 
0212  ' 
0213 
0214f- 
0215 
0216 
0217  -51’ 
0217.1 
0218 
0219  r 
022  0 
0221 
0222  e- 

0223  - 
0224 
0225  f 
0226  ■ 
0227 
0228 
0229  ' 
0230 
0231 
0232 
0233 
0234 
0235  " 
0236 
0237 
0238 
0239 
0240  V 
0241  j 
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GO  TO  6 

5  U  «  ABSF(Y/X) 

Q  »  (U-l. J/(U+1.) 

Q2  «  Q*0 

ALPHA  »  .78539815+Q*(  .99999612+Q2*(-.3;3»17376+Q2*(  .  1938784V+Q2* 
1  (-.13233510+Q2*(*07962632+02*(-'*03360&27-t-02**006Al2«lt’>)l)} 
ARCTAN  »  ABSF(R+ALPHA)  ” 

6  CONTINUE 
RETURN 
END 
DATA 


0242 

0243 

0244 


0245  0 

0246 

0247 

0248 

0249 

0250 


0251 

0252 


r»or»r>or>r»o  *  ♦ 


•  HERCO  PROJECT  P-106-9~  4  AHM  -  2-20-63  I TFAP-04 

XEQ 

FORTRAN 

INPUT  TIME  FUNCTION  APPROXIMATION  PROGRAM 
PROGRAMMED  FOR 

6593D  TEST  GROUP ( DEVELOPMENT )  -  AIR  FORCE  FLIGHT  TEST  CENTER  - 
EDWARDS  AIR  FORCE  BASE*  CALIFORNIA 
PROGRAMMED  BY 

HOUSTON  ENGINEERING  RESEARCH  CORPORATION 
UNDER 

CONTRACT  NUMBER  AF  04(611)  B199 
DIMENSION  IDRUN(2)  *ID0UT(  2)  •XCORDtSOOO)  vICTI  SOOD  )  *HEOIMG  (  36  )  * 

1  HW(5000) «QHW( 5000) tFWR ( 5000 ) »FW I ( 5000 ) »W( 5000 ) *XARRAY( 300) * 

2YARRAY(300)  .WARRAYOOO)  ’ 

EQUIVALENCE(HW*FWR) . (QHW^FWI ) .(XCORD*W) 

1  READ  INPUT  TAPE  5 . 100 .  (  HEDING  ( I  )  »  I «=1 . 36  ) 

100  FORMAT!  12A6/12A6/12A6) 

READ  INPUT  TAPE  5  *  101 . 1 DOUT ( 1 ) » IDOUT ( 2 ) •MBRPTO*BEFREQ»FIFR£Q* 
IFREQiN.SAHPLM' 

101  FORMAT(2I4.5X»I7.3F10.O.F15.8) 

READ  INPUT  TAPE  5 • 102 • ( W ( I ) *MW( 1 > •QHWI I > *  I -1 •NBRPTO ) 

102  FORMAT(10X.F10,0.2FI5.5) 

86  IBASE  =0 

23  IF! IBASE-NBRPT0)24.11.I1 

24  NPAGE  =  NBRPTO- IBASE 

Y  IF(NPAGE-200)25.26»26 

25  LIM2=NPAGE 

'  GO  TO  27 

26  LIM2  *  200 

*  27  CALL  HEDOUT! IBASE*HEOIMG) 

DO  199  I  =  1.LIM2 
K  =  I  +  IBASE 
WARRAYd  )  ■-  W(X) 

XARRAYII)  »  HW(X) 

199  YARRAYd)  *  QHWIK) 

WRITE  OUTPUT  TAPE  6 » 105 » I DOUT ( 1 ) » IDOUT ( 2 ) »  BEFREO*FI FREQ  . 

1FREQIN*!WARRAY(  I)  .XARRAY!  I)  .yARRAY!I)..I-l»LIM2) 

105  FORMATIIHO,  33HTRANSFER  FUNCTION  FOR  TRANSDUCER  *2 14 *6X»F10» 0. 

13H  -  iFlO.O.lOH  CPS  RANGE. 7X*Fi0,0»15H  CPS  INCREMENTS/ IH  / 
24(F8.0,4HCPS*»F7.3»3H  L'  .F8.3)) 

IBASE  =  IBASE  +  LIM2 
GO  TO  23 

11  READ  INPUT  TAPE  5  ♦  103 * lORUN (  1 )  DRUN ! 2  )  »NBRP TI  .PRE(3LM*OELTI  .DELTO* 
ITIME.CONAMP 

103  FORMAT(2I5*5X.I5.4F10*O.F1O.5) 

READ  INPUT  TAPE  5*  1 04 » ( XCQRD ( I ) »XT( I ) » I -I .NBRPT I ) 

104  FORMAT(14X»F6.0»F10.0) 

IFINBRPTI  -INBRPTI/2)»2)28.85#28 

L  85  NBRPTI  «  NBRPTI-1 
28  ISASE  «  0 


0001 
0002 
0003 
0004 
0005 
0006 
0007 
0008 
0009 
0010 
0011 
0012 
0013 
00.14 
0015 
0016 
00x7 
0018 
0019 
0020 
0021 
0022 
0023 
0024 
0025  . 
0026 
0027 
0028 
0029 
0030 
0031 
0032 
0033 
0034 
0036 
0036 
0037 
0038 
0039 
0040  . 
0041 
0042 
0043  , 
0044 
0045 
0046  . 
0047 
0048 
0049 
0058 


<  ) 


29  IF(  IBASE“NBRPTn31»35»35 

31  NPA6E  «  NBRPTI  -  IBASE 
IF(NPAGE  “  300)32*Jft33 

32  LIM2  »  NPAGE  ■ 

GO  TO  34 

33  LIM2  -  300 

34  CALL  HEDOUT  (ItHEOING) 

DO  168  I  -  It  LIH2 

K  ■=  I  +  IBASE 
XARRAYII)  4  XCORD(K) 

188  YARRAYd)  =  XT(K) 

WRITE  OUTPOT  TAPE  6 1 1 06 1 1 0RUN(  1 )  » IDROMJ 2  J.tNBRPTI  »D£LT 1 1 
1 ( XARRAY ( 1 1 » YARRA Y ( 1 ) i I >1 1 L IM2 ) 

106  FORMATdHO^ZSHTIME  RESPONSE  DATA  l.O*  «»21$tI12  rllH  POlMtS  AT  t 
lF10.3t  22(H  MICROSECOND  INTERVJ|iLS  /  6(2F10*0n 
IBASE  =  IBASE  +  LIM2  ’ 

GO  TO  29 

35  XCBASE  «  XCORD(l) 

XTBASE  «  XT( 1 ) 

DO  2  I  »  liNBRPTI 

XCORD(I)  <■  XCORDd)  -  XCBASE 

2  XT(I)  =  XTd)  -  XTBASE 
XTFIN  »  XTINBRPTI ) 

NBRW  =  (FREQLM  -  BEFREO ) /FREQIN 
IF(NBRW-(NBRW/2)»2)  At 3 #4 

3  NBRW  =  NBRW  +1 

4  FNBRW  -  NBRW 

FIFREQ  -  FRBRW*FREOIN 
INIT  -  50  ■  ■ 

FN  -  MBRPTI*NBRW 

ESTIME  »  FN*l*7E-5  +  TIME/DELT0*FNBRW»l*6E-5 
NBRT  »  TINE/DELTO»l.E6+  1. 

PRINT  50 1 1 1 DOUT ( 1 ) t IDOUT <  2  >  t lORUN ( 1 ) t IDRUNI 2 ) tEST IME 
501  FORMAT! IHOi  49HINPUT  TIME  FUNCTION  APPROXIMATION  FOR  TRANSDUCER  t 
12I4f  24H  AND  TIME  RESPONSE  I .O*  »2I5»/1H  t  38HMEXT  ON  LINE  PRINT  S 
IHOULD  OCCUR  WITHINiFlO.Of  8H  MlNUTES/lMl) 

FREQ  »  BEFREQ 
DO  49  J«ltNBRW 

ANGLE  *  6.283i893E-6*FREO*DELTI 
SUMOl  =0. 

S0MO2  -0. 

SUMEl  -  0. 

SUHE2’  •  0. 

LIM  -  N8RPTI  “  1 
D01331  ■  2tLIMt2 
FI  =  I-l 

TANGl  «  FI  *  ANGLE 

SUMOl  «  SUMOl  +  XT(1)*C0SF(TANG1) 

SUM02  -  SUM02  XT(  II«SINF(TAMG1 ) 

12  FI  »  I 


0051  0 

0052 

0053 

0054  0 

0055 

0056 


0057  0 

0058 

0059 


00600 

0061 

0062 


0063  0 
0064 
0065 
0066  0 
0067  ~ 
0068 


0069  0] 

0070 

0071 

0072  0 

0073 

0074 


0075 

0076 

0077 


0078  0 
0079 
0080 
0081  0 
0082 
0083 
0084  01 
0085 
0086 
0087  0 
0088 
0089 


0090  0 

0091 

0092 

0093  0 

0094 

0095 

0096  0 

0097 

0098 

0099  0 

0100 


© 
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TANGl  =  FI»ANGLE  0101 

SUMEi  =  SUMEl  +  XT( I+1)*C0SF(TANG1)  0102 

133  SUME2  =  SUME2  +  XT ( I+l ) *SINF ( TANGl )  0103 

SUMEl  =  SUMEl  -  XT(NBRPTI )*C0SF(TANG1 )  0104 

SUME2  =  SUME2  -  XT  (  NBRP T I  ) *S I NF  (  T  ANGl  )  0105  '' 

10  FN  =  NBRPTI  -  1  0106 

TANGl  =  FN*ANGLE  0107 

YOPYNl  =  XT(1)  #  XT(NBRPTI )*C0SF(TANG1 )  0108.' 

Y0PYN2  =  XT(NBRPTI  )*6IiMF(TAN&l)  0109 

SUMl  =  .33333333*(YOPYNl  +  4.*SUM01  +  2.*SUME1)  0110 

SUM2  =  .33333333*(Y0PYN2  +  4.*SUM02  +  2.*SUMt2)  0111: 

TANGl  =  ANGLE*! FN+.5)  0112 

FACT=1./SINF( .5*ANGLE) 

X  =(SUM1  -  .5*  XTFIN  *S I NF ( TANGl ) *F AC T ) *1 . E-6  0114, 

Y  =(-SUM2  -  .5*  XTFIN  *COSF (TANGl )*FACT  )*1 .E-6  011b 

AMP  =  SQRTF (X*X+Y*Y )  0116 

FMAG  =  AMP/HW(J)  0117,, 

FANG  =  ARCTAN(Y.X)  -  QHW ( J ) /57.2958  0118 

IF! FANG-6.2831853 ) 901 .901.900  0119 

900  FANG  =  FANG  -  6.2831853  0120* 

901  CONTINUE  0121  ' 

FWR!J)  =  FMAG*  COSF!FANG)  0122  ■' 

FWI!J)  =  FMAG*SINF(FAN:G)  0123~ 

QFANG  =  FANG*57.295a  0124" 

IF!INIT-50)  46,45.46  0125 

45  CALL  HEDOUT ! 1 .HEDING )  0126" 

WRITE  OUTPUT  TAPE  6 . 90  ,  I  DOUT  !  1  )  .  1 DOUT  !  2  )  .  I  DRUN  !  1 )  .  I  DKUN !  2  )  0127 

90  FORMATdH  ,  36HI NTERMEUI  ATE  RESULTS  FOR  TRANSDUCER  .214.  24H  AND  TI  0128 

IME  RESPONSE  I.D.  . 2 1 5 /8X , 4HFREQ . 9X  .AHMAGN . 8X  .  5HANGLE . 9X . 4HREAL .  0129'' 

29X.4H.IMAG.  12X  .IHX  .  12X.  IHY  .9X  .4HSUivll.9X.4HSUM2/lH  )  0130 

INIT  =  0  0131 

46  WRITE  OUTPUT  TAPE  6 . 8u . FREQ . FMAG  .QFANG .F WR ! J )  . FW I ! J ) .X . Y .SUM  1 . SUM2  0132  ’' 

80  FORMATdH  .  F 1 1 . 0 , 8  El  3. 4  )  0133 

^NIT=INIT+1  0134 

49  FREQ  =  FREQ  +  FREQIN  0135" 

INPUT  TIME  FUNCTION  COMPUTATION  0136 

INIT=50  0137 

T  =  0.  0138'' 

601  SUMOl=0.  0139 

SUM02  =  0.  m  -  0140 

SUME1  =  0.  0141 

SUME2=0.  0142' 

W=BEFREQ*6. 2831853  0143  i 

LIM  =  NBRW  -  1  ^  0144  '  ! 

DO  602  1=2. LIM. 2  0145 

FI1=I-1  0146 

FI2=I  0147 

ARG1=!FREQIN*6.2831853*FI1+W)*T  0148 

ARG2=!FREQIN*6.2831853*FI2+W)*T  0149 

SUM01  =  SUM01+FWR! I )*COSF( ARGl )  0150  “ 


SUME1=SUME 1+FWR ( I +1 ) *COSF ( ARG2 ) 

SUM02-‘:0'MG.;+|-WI  ( I  )*SINF(ARG1  ) 

602  SUME2=SUME2+FWI ( I +1 ) *SI NF ( ARG2 ) 

T1=BEFREQ*6.2831853*T 
T2=FIFREQ*6.2831853*T 
T3=FWR(NBRW  )*COSF(T2) 

T4=FWI(NBRW  )*SINF(T2) 

SUME1=SUME1-T3 
SUME2=SUME2-T4 
Y0PYN1=FWR{1)*C0SF(T1 )+T3 

YOPYN2=FWI ( 1 )*SINF( T1 )+T4  • 

DEL=FREQIN*( 1.3333333  )#DELTI 
SUM1=DEL*( Y0PYN1+4.#SUM01+2.*SUME1 ) 

SUM2=-DEL*( Y0PYN2+4.*SUM02+2.*SUME2 ) 

IF( INIT-50)606.603.6U3 

603  CALL  HEDOUK l.HEDING) 

WRITE  OUTPUT  TAPE  6 . 605 . I DOUT ( 1 ) . I DOUT ( 2 ) ♦ I DRUN ( 1 ) . I DRUN ( 2 ) . FREU Lm 

1  .NBRT.DELTO.CONAMP.SAMPLM 

605  FORMAT! lHUii6X.44HC0MPUTED  INPUT  TIME  FUNCTION  FOR  TRANSDUCER  .214 
1.24H  AND  TIME  RESPONSE  I.D.  .2I5/1H  »19HCUT~OFF  FREQUENCY  =.F10.0. 

2  4H  CPS.  I1U.17H  OUTPUT  POINTS  AT.F10.3  .22H  MICROSECOND  INTERVAL 
3S/20H  CONVERSION  FACTOR  = . F 1 0 . 5 . 1  OX . 39HALLOWABLE  APPROXIMATION  TIM 
4E  INTERVAL  =.F10.8.8H  SECONDS/ 

513X.17HTIME-MILLISECONDS. 15X . 15HCOSI NE  I  NTE6RAL. 17X. 13HSINE-INTEGR 
6AL.17X.13HTIME  RESPONSE/IH  ) 

INIT=0 

606  TYME=T*1.E3 

K  =  T/DELTI*1.E6  +  1. 

SUMl  =  SUM1*C0NAMP 
SUM2  =  SUM2*CONAMP 
TIRESP  =  XT(K)*CONAMP 

WRITE  OUTPUT  TAPE  6 . 7 02 . T YME .SUMl »SUM2 .T IRESP 
702  FORMAT(F30.3.3F30.5 ) 

INIT=INIT+1  : 

IF(T-TIME)610.611.611 

610  T=T+DELTO*l,E-6 
GO  TO  601 

611  PRINT  600.IDOUT(1) .IDOUT(2) .lURUN(l) .IDRUN(2  ) 

600  FORMAT! IH  .  49HINPUT  TIME  FUNCTION  APPROXIMATION  FOR  TRANSDUCER  » 
1214.  19H  and  TIME  RESPONSE  .2I5.10H  COMPLETED/ IHl ) 

GO  TO  1 
END 

FORTRAN 

SUBROUTINE  HEDOUT !  M.HEDING) 

DIMENSION  HEDING!36) 

IF!  M  )  802.801.802 

801  IPAGE  =  M 

802  IPAGE  =  IPAGE  +  1 

WRITE  OUTPUT  TAPE  6.  8 03 . 1  PAGE . I HEDI NG ( I ) . I  =  1 .36 ) 

803  F0RMAT(1H1.57X»5HPAGE  . I  2 » /24X . 12A6/24X. 12 A6 /24X . 1 2 A6 ) 


0151  , 
0152  ’ 
0153 
0154 
0155-,^' 
0156 
0157 
oiss;-; 
0159 
0160 
0161  ' 
0162 
0163 
0164;'*' 
0165 
0166 
01670 
00168 
0169 
0170'“'; 
0171 
0172 
0173  ~ 
0174 
0175 
0176  . 
0177 
0178 
0179  ‘ 
0180  " 
0181 
0182  ■ 
0183 
0184 
0185  ■*'■ 
0186 
0187 
0188 
0189 
0190 
0191 
0192 
0192.1 
0193~': 
0194 
0195 
0196  ■' 
0197  ' 

0198 
0199' 
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RETURN 
t  END 

*  FORTRAN 

FUNCTION  ARCTAN(Y*X) 

IF  (X)  1.2,3 
*  1  IF(Y)11,12,13 

11  R=  3.1415927 
GO  TO  5 

12  R  -  3.1415927 

^  GO  TO  4 

13  R  *  -3.1415927 
GO  TO  5 

2  IFCY)  21.22t23 

21  R  =  4.7123»e9 
GO  TO  4 

22  R  »  0. 

60  TO  4 

23  R  =■  1.5707963 
GO  TO  4 

3  IF  (Y)  31.32,33 

31  R  =  -6;  28314154 

GO  TO  5  ' 

32  R  =  0. 

GO  TO  4  • 

33  R  =  0. 

'*  GO  TO  5 

4  ARCTAN  *  R 

'  GO  TO  6 

5  U  *  ABSF<Y/X) 

•  Q  «=  (U-1.  ) /(U+1.  ) 

02  =  Q»Q 

ALPHA  =  .78539815+Q*( .99999612+Q2*(-.33317376+02*( .19e07e69+Q2» 
1 ( - . 1 32  3  3  5 1 0  +Q2 * ( . 0  79626 3 2  +Q2  * ( - . Q 3 360627 +02 * . 0 068 1 24 1 ) »! ) ) ) 
ARCTAN  =  ABSF(R+ALPHA) 

6  CONTINUE 
RETURN 
END 
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